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Summary 

A previously developed analytical fonnulation has been modified in order to more accurately 
account for the effects of hydrostatic stresses on the nonlinear, strain rate dependent deformation 
of polymer matrix composites. State variable constitutive equations originally developed for 
metals have been modified in order to model the nonlinear, strain rate dependent defonnation of 
polymeric materials. To account for the effects of hydrostatic stresses, which are significant in 
polymers, the classical h plasticity theory definitions of effective stress and effective inelastic 
strain, along with the equations used to compute the components of the inelastic strain rate 
tensor, are appropriately modified. To verify the revised formulation, the shear and tensile 
deformation of two representative polymers are computed across a wide range of strain rates. 
Results computed using the developed constitutive equations correlate well with experimental 
data. The polymer constitutive equations are implemented within a strength of materials based 
micromechanics method to predict the nonlinear, strain rate dependent defonnation of polymer 
matrix composites. The composite mechanics are verified by analyzing the deformation of a 
representative polymer matrix composite for several fiber orientation angles across a variety of 
strain rates. The computed values compare well to experimentally obtained results. 

Introduction 

NASA Glenn Research Center has an ongoing research project to develop computational 
methods for the analysis of polymer matrix composites subject to high strain rate impact loads. 
Under these types of loading conditions, the deformation of the composite can be highly strain 
rate dependent and nonlinear, which must be accounted for within the analytical model. The 
deformation of polymer composites is often assumed to be linear elastic and independent of 
strain rate in transient dynamic finite element codes used for impact analysis (ref. 1). 

Polymers are known to have a strain rate dependent deformation response that is nonlinear 
above about one or two percent strain. Traditionally, viscoelasticity models have been used to 
capture this behavior (ref. 2). However, there has been an interest in the research community in 
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using constitutive equations developed for metals, based on plasticity and viscoplasticity theory, 
to model the nonlinear, strain rate dependent deformation of polymers and polymer matrix 
composites. For example, Sun and co-workers (refs. 3, 4) developed a macromechanical, 
transversely isotropic plasticity theory to analyze the nonlinear deformation of polymer 
composites. Bordonaro (ref. 5) adapted the viscoplasticity theory based on overstress, originally 
developed for metals, to analyze the nonlinear deformation of Nylon 66. Pan and co-workers 
(refs. 6, 7) and Hsu, et al (ref. 8) developed viscoplasticity theories to analyze the nonlinear 
deformation of polymers including the effect of hydrostatic stresses. In these studies, the effects 
of hydrostatic stresses on the nonlinear deformation, which unlike in metals are significant for 
polymers (ref. 9), were accounted for by applying variations of the Drucker-Prager yield 
criterion (ref. 10) to modify the definitions of the effective stress and the effective inelastic 
strain. A preliminary effort to account for the effects of hydrostatic stresses on the high strain 
rate defonnation of polymers has been carried out by the authors of this report (ref. 11). 

In this paper, constitutive equations with an associated flow rule are developed to analyze the 
nonlinear, strain rate dependent deformation of polymer matrix materials in which the effects of 
hydrostatic stresses on the nonlinear deformation response are properly accounted for. The 
equations developed in previous work (ref. 11) were based on a non-associated flow rule and 
were found to have difficulties in accurately predicting the multiaxial stress states found in the 
matrix of a composite material. To calibrate and check the accuracy of the revised equations, 
two representative toughened epoxy polymers (PR520 and 977-2) are characterized, and the 
tensile and shear deformation of these materials are modeled for a wide range of strain rates. 

The polymer constitutive equations are implemented within a mechanics of materials based 
micromechanics method incorporating a fiber substructuring approach. This capability pennits 
the analysis of the strain rate dependent, nonlinear deformation of polymer matrix composites. 
By using this approach, first the matrix constituent is characterized, and its deformation response 
is computed using the polymer constitutive equations, and then micromechanics equations are 
applied along with the fiber properties to compute the effective properties and deformation 
response of the composite. A justification for applying this type of micromechanics approach is 
that conducting tensile and shear tests on the bulk polymer is much easier than conducting 
equivalent tests on the composite, particularly at high strain rates and for composites with 
complicated fiber layups. Therefore, developing a methodology which reduces the amount of 
testing that is required on the composite, but still can yield strain rate dependent, nonlinear 
deformation data which can be used by a finite element code, is desirable. In this methodology, 
the in-situ matrix properties and defonnation response are assumed to be equivalent to those of 
the bulk polymer. This assumption is useful in that it is simpler, due to the number of material 
constants required for the polymer equations and the procedures used to obtain the constants, to 
conduct experiments on the pure resin and to determine the resin material constants with this 
data than it is to try to back out the resin properties from composite test data. However, if in 
comparing test data obtained from composite specimens to analytical predictions it appears that 
the bulk matrix properties do not accurately reflect the in-situ state of the matrix, the polymer 
properties can always be appropriately adjusted. To verify the micromechanics methodology, 
the tensile deformation of a representative carbon fiber reinforced composite, with one of the 
polymers used to correlate the constitutive equations (977-2) as its matrix constituent, is 
analyzed for several fiber orientation angles across a wide range of strain rates. The analytical 
predictions are compared to experimentally obtained results. 


NASA/TM— 2003-212382 


2 



Polymer Constitutive Equations 

Overview 

To analyze the nonlinear, strain rate dependent deformation of the polymer matrix materials, 
the Bodner-Partom viscoplastic state variable model (ref. 12), which was originally developed to 
analyze the viscoplastic deformation of metals above one-half of the melting temperature, is 
modified. In state variable models, a single unified strain variable is defined to represent all 
inelastic strains (ref. 13). Furthermore, in the state variable approach there is no defined yield 
stress. Inelastic strains are assumed to be present at all values of stress, the inelastic strains are 
just assumed to be very small in the “elastic” range of deformation. State variables, which 
evolve with stress and inelastic strain, are defined to represent the average effects of the 
deformation mechanisms. 

For this study, temperature and moisture effects have been neglected, as only room 
temperature data are currently available. These effects will be examined in future efforts. The 
nonlinear strain recovery observed in polymers on unloading is not currently simulated, and 
small strain theory is assumed to apply. Phenomena such as creep, relaxation and high cycle 
fatigue are not currently accounted for within the equations. However, since hydrostatic stress 
effects have been found to be important in the analysis of high cycle fatigue (ref. 14), correctly 
accounting for the effects of hydrostatic stresses will be critical if in the future fatigue analysis 
capabilities are added to the polymer constitutive equations. All of the nonlinearity and strain 
rate dependence observed in the polymer is assumed in this work to be due to inelastic 
deformation, where in reality the nonlinearity could be due to a mixture of deformation and 
damage. 

An important point to note for this work is that in the original Bodner model as applied to 
metals, each of the state variables and material constants had a fairly explicit link to specific 
deformation mechanisms. In the application of the equations to the analysis of polymers, the 
relationship of the various constants and variables to specific deformation mechanisms is not as 
well defined. An internal stress state variable is defined to represent phenomelogically the 
resistance to molecular flow that takes place during the polymer defonnation. An additional 
state variable is defined which represents phenomelogically the level of influence of the 
hydrostatic stresses on the inelastic defonnation. The remainder of the material constants are 
characterized in order to match the experimental deformation response of the polymer, and are 
not tightly linked to specific defonnation mechanisms. However, the equations developed here 
employ a fairly simple formulation with a minimum of easily obtained material constants, and as 
will be discussed later appear to simulate the nonlinear, strain rate dependent deformation of the 
polymer fairly accurately both qualitatively and quantitatively. 

Associated Flow Equation and Evolution Equations 

In order to derive the associated flow equation for the components of the inelastic strain rate 
tensor for polymers with hydrostatic stress effects included, a procedure similar to that employed 
to derive the Prandtl-Reuss equations in classical plasticity is utilized (refs. 10, 13). The 
formulation employed in this study is based on that used by Pan and co-workers (refs. 6, 7) and 
Hsu, et al (ref. 8). First, an inelastic potential function based on the Drucker-Prager yield 
criterion (ref. 10) is defined as 


I ^[^2 ®®kk 


( 1 ) 
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where J 2 is the second invariant of the deviatoric stress tensor. The variable Okk is the sum of the 
normal stress components (assuming the usual rules of indicial notation (ref. 13)) and is equal to 
three times the hydrostatic stress. The variable a is a state variable which controls the level of 
the hydrostatic stress effects. The tenn “aakk” incorporates the effects of hydrostatic stresses 
into the inelastic potential function. The components of the inelastic strain rate tensor, , are 

assumed to be proportional to the derivative of the potential function with respect to the 
components of the stress tensor, Oy, as follows 



JL 

da IJ 


( 2 ) 


where A is a scalar rate variable. By taking the specified derivative of the potential function, the 
result becomes 


3/ _ 
u 2^/77 


+ aS tj 


(3) 


where Sy are the components of the deviatoric stress tensor and 8y is the Kronecker delta. 

To solve for A , Equation 3 is substituted into Equation 2, and by multiplying the resulting 
equation by itself, the following expression is obtained 


E 1 E 1 


4J, 


0 C 

+ ^UL a + a 2 8 li 8 lj 


U 2 


(4) 


Equation 4 can be simplified by applying the following easily proven identities 


S,Sy = 2 J 2 

Sy8y = S, = 0 

M = s« = 3 


which leads to the following expression 


s'jE^X- 


v47 


+ 3 a" 


=x- 


f l + 6er 'l 


which can be solved for A to obtain 


(5) 


( 6 ) 


A = 


V 


2 E 1 E 1 


Vl + 6cr 2 


( 7 ) 
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By substituting Equation 3 and Equation 7 into Equation 2, the associated flow rule can be 
specified 


4 = 






Vl + 6a 2 


S ij 


2^ 


( 8 ) 


The next step in defining the flow rule is to define the effective stress and the effective 
inelastic strain rate. The effective stress, c e , is defined as follows 

o e = Sf = + V3 acj kk . (9) 

Under pure shear loading, or for cases where the state variable a equals zero, the hydrostatic 
stress tenn becomes zero and the effective stress definition reduces to the classical definition 

of y/3J 2 . To detennine the effective inelastic strain rate, e { , the principal of the equivalence of 
the inelastic work rate is employed. With the aid of Equations 2 and 3, the inelastic work 
rate, W 1 , can be expressed by the following 


W 1 








2^ +aSiJ 


U = 


( 10 ) 


After applying the effective stress definition given in Equation 9 and simplifying, it can be 
shown that 


A = V3. 


el 


( 11 ) 


which by comparison to Equation 7 leads to a definition for the effective inelastic strain rate 


el = 


1 


2 ./ .7 

— ff 

3 lJ lJ 


yll + 6a 2 


( 12 ) 


By substituting Equation 1 1 and Equation 3 into Equation 2, an alternative definition for the 
components of the inelastic strain rate tensor can be given as 


e‘j = V3^ 7 


W, +aS “ 


(13) 


By utilizing Equation 13, it can be shown that an alternative definition of the effective inelastic 
strain rate can be given as 
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( 14 ) 


* =ei =J 


e 1 . = e I . -e 1 8 

ij ij m ij 


where e' e is the effective deviatoric inelastic strain rate and e’ m is the mean inelastic strain rate, 
which matches the effective inelastic strain rate definition given by Pan and co-workers (refs. 6, 
7). 

By following a similar format to the Bodner-Partom model (refs. 12, 13), the following 
definition is specified 



D 0 exp 


/ \ 2 n 

r z A 


\°ej 


(15) 


where D 0 and n are material constants. D 0 represents the maximum inelastic strain rate, and n 
controls the rate dependence of the material. Z is an isotropic state variable which represents the 
resistance to molecular flow (internal stress). By substituting Equation 15 into Equation 13, the 
final form of the flow equation is determined to be 


K = 2 D a exp 


/ \ 2 n 

r z A 




W> +aS ’ 


(16) 


Note that the elastic components of the strain rate must be added to the inelastic strain rate to 
obtain the total strain rate. 

The rate of evolution of the internal stress state variable Z and the hydrostatic stress effect 
state variable a are defined by the equations 


Z=g(Z l -Z)ei 

(17) 

a = q(a l -a)e , e 

(18) 


where q is a material constant representing the “hardening” rate, and Zi and oq are material 
constants representing the maximum values of Z and a, respectively. The initial values of Z and 
a are defined by the material constants Z Q and a 0 . The tenn e[ in Equations 17 and 18 

represents the effective deviatoric inelastic strain rate, which was defined in Equation 14. An 
important point to note is that in the original Bodner model (ref. 12), the inelastic work rate was 
used instead of the effective inelastic strain rate in the evolution equation for the internal stress 
state variable. However, for this work the inelastic strain rate was deemed easier to work with 
from both computational and characterization points of view, particularly in the incorporation of 
hydrostatic stress effects. Furthermore, the effective inelastic strain rate has been used in other 
state variable constitutive models (ref. 13). Since hydrostatic stress effects were not considered 
in the original Bodner model (ref. 12), the evolution equation for a is new to this work. The 
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state variable a is assumed to evolve in the same manner as the state variable Z. As will be 
pointed out later, by using this assumption the value of q used in Equation 17 will be the same as 
the value of q used in Equation 18. 

Detennination of Material Constants 

The material constants that need to be determined include Do, n, Z 0 , Zi, a 0 , oq, and q. The 
procedure to be used is summarized here. More details on the general approach can be found in 
Stouffer and Dame (ref. 13) and Bodner (ref. 12). The values of D 0 , n, and Zi are characterized 
as follows using Equation 16. The value of D 0 is currently assumed to be equal to a value of 
10 4 times the maximum applied strain rate, which correlates with the maximum inelastic strain 
rate. Equation 16 is simplified to the case of pure shear loading, leading to the following 
expression 


f 

2 


D 0 exp 



r 

Id 


(19) 


where y l is the engineering shear strain rate, x is the shear stress, and the remainder of the tenns 
are as defined earlier. The case of pure shear loading is used to characterize the equations for 
two reasons. First, hydrostatic stress effects are not present in the case of pure shear loading, so 
those material constants do not enter into the picture, simplifying the initial characterization. 
Second, since polymers tend to be more ductile in shear than in uniaxial tension, the shear stress- 
shear strain curves obtained experimentally are more likely to display a defined “saturation” 
stress, which as shown below is crucial for determining the material constants. The “saturation”, 
or “yield”, stress is defined as the stress level where the inelastic strain rate equals the total 
applied strain rate. This condition occurs at the point where the stress-strain curve flattens out 
and becomes horizontal. 

Next, Equation 19 is rearranged as follows 


- 2 In 


(f ) 


( z 

v2A> J 




2 n 


( 20 ) 


and the natural logarithm of both sides of the resulting expression is taken. The values of the 
inelastic shear strain rate, shear stress, and state variable Z at “saturation” are substituted into the 
resulting expression, which yields the following 


In 

-2 In 

r To ) 



V 2 A) 


2 n ln(Zj ) - 2 n ln( V3r s ) 


( 21 ) 


where x s equals the saturation shear stress, f o is the constant applied total engineering shear 
strain rate in a constant strain rate shear test, and the remaining terms are as defined earlier. 

The required constants are determined from a set of shear stress-shear strain curves obtained 
from constant strain rate tests. Each curve in this set is obtained at a different constant strain 
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rate. Data pairs of the total strain rate and saturation shear stress values from each curve are 
taken. For each strain rate, the data values are substituted into Equation 21, and represent a point 
on a master curve. The number of points in the master curve equals the number of strain rates at 
which shear tests were conducted. A least squares regression analysis is then performed on the 
master curve. As suggested by Equation 21, the slope of the best- fit line is equal to -2n. The 
intercept of the best-fit line is equal to 2n(ln (Zi)). 

To detennine the value of Z 0 , first Equation 20 is rearranged as follows 




( • I \ 

\ 

z = 

-2 In 

f 7 1 




2 D 1 



V 

V o J 

J 


& i 


( 22 ) 


where all of tenns are as defined earlier. To find the value of Z 0 from Equation 22, the value of 
the shear stress where the stress-strain curve becomes nonlinear for a particular constant strain 
rate shear test is used for the value of x. The point where the stress-strain curve becomes 
nonlinear is defined as the approximate point where the curve appreciably deviates from a linear 
extrapolation of the initial data. The value of y l is set equal to the approximate inelastic shear 
strain rate when the stress-strain curve becomes nonlinear. The shear strain rate used in the 
constant strain rate test divided by 100 was found by trial and error to approximate this value 
reasonably well. Using this data, Equation 22 is solved for Z, which is assumed to be equal to 
the value of Z 0 . Using the data from the lowest strain rate test available has been found to give 
adequate values of Z 0 . However, the calculations can be made using data from all the available 
strain rates, and an average taken if required to obtain the value of the constant. 

To detennine the value for q for Equations 17 and 18, first Equation 17 is integrated for the 
case of pure shear loading, resulting in the following relation 


Z = Z 1 -(Z l -Zjexp 


f -q / 

s 7 ) 


(23) 


where y 1 is the inelastic shear strain. At saturation, the value of the internal stress Z is assumed 
to approach Z\, resulting in the exponential tenn approaching zero. Assuming that saturation 
occurs when the following condition is satisfied 


exp 


r -q 

S 7s 


= 0.01 


(24) 


the equation is solved for q, where y s J is the inelastic shear strain at saturation. If the inelastic 
shear strain at saturation is found to vary with strain rate, the parameter q is computed at each 
strain rate and regression techniques are utilized to detennine an expression for the variation of 
q. If Equation 18 is integrated, an expression similar to Equation 23 is obtained. At saturation, 
the value of a is assumed to approach a i, so an equation identical to Equation 24 is obtained, 
which would lead to the same value for q. Therefore, identical values of q are used in Equation 
17 and Equation 18. 


NASA/TM— 2003-212382 



To obtain the values of oq and a 0 , Equation 9 is used in combination with stress-strain data 
from constant strain rate uniaxial tensile tests and constant strain rate shear tests. The primary 
assumption used at this point (and assumed implicitly in Equation 9) is that the effective stress at 
saturation under uniaxial tensile loading at a particular strain rate is equal to the effective stress 
at saturation under pure shear loading at the same equivalent strain rate. Likewise, the effective 
stress at the point the stress-strain curve becomes nonlinear under tensile loading is equal to the 
effective stress at the point the stress-strain curve becomes nonlinear under shear loading. 
Therefore, assuming the value of a at saturation is equal to oq, and the value of a at the point the 
stress-strain curve becomes nonlinear is equal to a Q , the following equations are obtained for the 
case of having data from uniaxial tension tests and pure shear tests 


7 s (l + V3a,)=V3r s 

(25) 

„,(l + V3a 0 )=V3r (I , 

(26) 


where o s and x s are the tensile and shear stresses at saturation, respectively, and c n i and T n i are 
the tensile and shear stresses at the point where the respective stress-strain curves become 
nonlinear. The required constants can then be determined from these equations. The values of 
the material constants are assumed to be rate independent, so the results from only one strain rate 
need to be used to find the needed parameters. In practical application of the methodology, the 
uniaxial tension and pure shear tests used do not have to be at the exact same effective strain 
rate. As long as the effective strain rates from the two tests are approximately equal, the values 
obtained have been found to be valid. 

Simulation of Strain Rate Dependent Polymer Deformation 

To demonstrate the ability of the developed constitutive equations to correctly simulate the 
hydrostatic stress dependent deformation response of polymers, two representative toughened 
epoxies, PR520 and 977-2, were analyzed. Details of the experimental procedure and 
experimental results will be given in a future report. However, a summary is given here for 
completeness. Longitudinal tensile tests and pure shear tests were conducted at room 
temperature on the materials at strain rates of about 5xl0' 5 /sec, 1 /sec and 400 /sec. The low and 
moderate strain rate tests were conducted using an Instron hydraulic testing machine. The high 
strain rate tests were conducted using a split Hopkinson bar. Engineering stress and engineering 
strain were measured until failure. 

Shear stress-shear strain curves for PR520 obtained under pure shear loading are shown in 
Figure 1 for each of the strain rates examined, while tensile stress-strain curves are shown in 
Figure 2. Shear stress-shear strain curves for 977-2 are shown in Figure 3, and tensile stress- 
strain curves are shown in Figure 4. Both materials exhibit a strain rate dependent, nonlinear 
deformation response under both shear and tensile loading. At high strain rates, the sharp 
increase in stress at the beginning of the loading with negligible increase in strain observed for 
both materials under shear loading and for PR520 under tensile loading is most likely the result 
of a lack of stress equilibrium at the start of loading. Due to this characteristic of the high strain 
rate shear curves, the initial elastic moduli at high strain rates for both materials were determined 
by approximations and extrapolations of the data. The oscillations seen in the tensile response of 
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977-2 at high strain rates are most likely due to the specimen geometry, leading to stress waves 
being visible in the response. The tensile specimen design for high strain rate testing was 
changed between the times when the 977-2 tensile data were obtained and when the PR520 data 
were obtained, which is why the oscillations do not appear in the PR520 data. The failure 
stresses under tensile loading for 977-2 appear not to vary with strain rate, and the failure stress 
for PR520 does not appear to increase in going from the low to the medium strain rates. 
Preliminary investigations indicate that the measured failure stresses at the higher strain rates 
may be artificially low due to the presence of strain gages on the specimen and the geometry of 
the specimen. Further details of these studies will be given in a future report. The material 
constants for both polymers were detennined using the procedures described earlier in this report 
and are listed in Table 1. 

The shear stress-shear strain curves computed for all three strain rates, along with the 
experimental results for comparison, are shown in Figure 1 for PR520 and Figure 3 for 977-2. 
Overall, the computed results correlate well with the experimental values for all strain rates for 
both materials. Specifically, the nonlinearity and rate dependence of the experimental results are 
captured qualitatively, and the quantitative match between the experimental and computed 
results is reasonably good. The high strain rate results are somewhat under predicted for both 
materials (particularly for PR520) at lower strains, but this is most likely due to the fact that in 
the experiments the initial stresses increased significantly with a negligible increase in strain, and 
thus the initial modulus of the material was computed by extrapolating from data obtained after 
the strain became non-negligible. 

The tensile stress-strain curves computed for all three strain rates, along with the 
experimental results for comparison, are shown in Figure 2 for PR520 and Figure 4 for 977-2. 

For both materials, qualitatively the nonlinearity and rate dependence of the experimental results 
is captured. Quantitatively, for PR520 at the medium strain rate the stresses in the nonlinear 
range are somewhat under predicted and for 977-2 at the low strain rate the stresses in the 
nonlinear range are somewhat over predicted, but otherwise there is a good correlation between 
the experimental and computed curves. There is not currently a good explanation for the cause 
of the discrepancies, but the maximum total error is still less than ten percent. The high strain 
rate tensile stresses for PR520 are somewhat under predicted at lower strains, but this is most 
likely due to the methods used to detennine the initial elastic modulus and the fact that the 
experimental results showed a sharp increase in stress with negligible increase in strain at low 
strains. The high strain rate tensile results computed for 977-2 bisect the oscillating 
experimental curve reasonably well. For the tensile results, the important point to note is that the 
material constants were primarily computed using the shear data, and the comparison of the 
computed results to the experimental tensile data is overall still reasonably good. 

To further explore the significance of properly accounting for the hydrostatic stress effects in 
the analysis, the tensile stress-strain curve for PR520 at the low strain rate of 5x1 O' 5 /sec is once 
again considered. In Figure 5, the experimental stress-strain curve, along with the original 
computed curve (“Alpha Included” in the figure), is presented. A tensile curve computed in 
which the hydrostatic stress effects are neglected (a o =ai=0, “Alpha=0” in the figure) is also 
shown. As can be seen in the figure, when the hydrostatic stress effects are neglected in 
computing the tensile response, the stresses are significantly over predicted as compared to the 
experimental results, which indicates that the hydrostatic stress effects are significant for 
polymers and must be accounted for within an analysis. 
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Composite Micromechanical Modeling 

Micromechanics techniques are used to predict the effective properties and deformation 
response of the individual plies in a composite laminate. The effective properties and 
deformation response are computed based on the properties of the individual constituents. 
Lamination theory can then be used to compute the effective deformation response of the entire 
composite. The constitutive equations described above have been implemented within a 
micromechanics method in order to enable the prediction of the nonlinear, strain rate dependent 
deformation response of polymer matrix composites with the effects of hydrostatic stresses 
incorporated into the analysis. The micromechanics method has been described extensively in 
Goldberg (ref. 15). A summary of the methodology will be given here. 

For this work the unit cell, the smallest material unit for which the response can be 
considered to be representative of the entire composite ply, is defined to consist of a single fiber 
and its surrounding matrix. Due to symmetry, only one-quarter of the unit cell is analyzed. The 
composites are assumed to have a periodic, square fiber packing with a perfect interfacial bond. 
The fibers are assumed to be transversely isotropic and linear elastic with a circular cross- 
section. The matrix is assumed to be isotropic, with a rate dependent, nonlinear defonnation 
response computed using the equations described in the previous section of this report. As 
mentioned previously, the in-situ matrix properties and defonnation response are assumed to be 
equivalent to those of the bulk polymer. 

The unit cell is divided up into an arbitrary number of rectangular, horizontal slices of equal 
thickness, as is shown in Figure 6. Similar approaches have been used by researchers such as 
Whitney (ref. 16), Greszczuk (ref. 17) and Mital, et al (ref. 18). Each slice is assumed to be in a 
state of plane stress and classical laminate theory is assumed to be applicable at this scale. The 
top and bottom slices in the unit cell are composed of pure matrix. The remaining slices are 
composed of two subslices; one composed of fiber material and the other composed of matrix 
material. For the slices containing both fiber and matrix, the out-of-plane stresses can be 
nonzero in individual subslices, but the volume average of the out-of-plane stresses must be 
equal to zero. By using this approach, the behavior of each slice is decoupled, and the response 
of each slice can be determined independently, which significantly reduces the level of 
complexity in the analysis. Laminate theory is used to obtain the effective response of the unit 
cell. In the future, the ability to incorporate transverse shear stresses and transverse pressures 
will be added to the model as these loading conditions might be significant in an impact analysis. 
In this case, the out-of-plane stresses would be constant throughout the unit cell. However, if the 
capability for out-of-plane loading is added to the model, adjustments to the model might need to 
be made to account for the fact that by dividing the unit cell using horizontal slicing, the overall 
unit cell properties might be orthotropic instead of transversely isotropic. In reality, the unit cell 
properties should be transversely isotropic, but since in the current study only in-plane stresses 
and loadings are assumed to be present, and the unit cell is assumed to be in a state of plane 
stress, this possible deficiency in the current model should not be significant. 

The thickness, fiber volume ratio and thickness ratio (the ratio of the slice thickness to the 
total unit cell thickness) for each slice can be determined using the composite fiber volume ratio 
and geometric principles. The unit cell is assumed to measure one unit in length by one unit in 
height. The first step is to compute the area of the cross-section of the fiber within each slice. 
The overall diameter of the fiber (df) is related to the fiber volume fraction of the overall 
composite (Vf) through the following relationship 
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and this tenn can be used along with the standard geometric definition of the radius of a circle to 
compute the horizontal coordinate of any point on the outer surface of the fiber in tenns of the 
liber volume fraction and the vertical coordinate. The area of the portion of the liber contained 
within each slice (A/) of the one-quarter of the unit cell which is analyzed can be computed by 
integrating the resulting expression between the vertical (z) coordinates of the top and bottom of 
slice “i” 



which is also the equivalent area of the rectangular liber slice in the one-quarter of the unit cell 
which is analyzed. The fiber diameter and liber area within each slice can be used to compute 
the fiber volume fraction and thickness ratio of each slice. 

The effective properties and effective inelastic strains of each slice are computed 
independently. The responses of each slice are combined using laminate theory to obtain the 
effective response of the corresponding lamina. Micromechanics equations are developed for 
those slices composed of both liber and matrix material. The stresses in the slices composed of 
pure matrix can be computed using the matrix elastic properties and inelastic constitutive 
equations. The standard transversely isotropic compliance matrix (or isotropic in the case of the 
matrix) is used to relate the local strains to the local stresses in the fiber and matrix. Each slice is 
assumed to be in a state of plane stress on the global level, but out-of-plane normal stresses can 
exist in each subslice. 

Along the liber direction, the strains are assumed to be uniform in each subslice, and the 
stresses are combined using volume averaging. The in-plane transverse nonnal stresses and the 
in-plane shear stresses are assumed to be uniform in each subslice, and the strains are combined 
using volume averaging. The out-of-plane strains are assumed to be uniform in each subslice. 
The volume average of the out-of-plane stresses in each subslice is assumed to be equal to zero, 
enforcing a plane stress condition on the global level for the slice. 

An orthotropic compliance matrix is used to relate the strains to the stresses in each 
constituent. The inelastic strain components are added to the standard orthotropic elastic 
compliance matrix to incorporate inelasticity into the constitutive relations. For the fiber, which 
is assumed to be linear elastic, these components are neglected. 

By combining the unifonn stress and uniform strain assumptions with the constituent stress- 
strain relations, a system of four simultaneous equations results that can be solved for the 
unknown stresses in the subslices. The total strains and subslice inelastic strains are considered 
to be the known values in solving this problem. By substituting the subslice stresses back into 
the equations defining the uniform stress assumptions, the effective elastic constants and 
effective inelastic strains can be computed. By applying classical laminate theory at this point, 
the effective stiffness matrix and effective inelastic strains for the unit cell are computed. 
Laminate theory is applied once again to obtain the effective properties and force resultants due 
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to inelastic strains for the multilayered composite laminate. Further infonnation on all of these 
procedures can be found in Goldberg (ref. 15). 

Simulation of Strain Rate Dependent Composite Deformation 

To verify the micromechanics equations, and the implementation of the polymer constitutive 
equations within them, a series of analyses have been carried out using a representative polymer 
matrix composite system that exhibits a strain rate dependent, nonlinear deformation response. 
The material examined consists of carbon IM7 fibers in the 977-2 toughened epoxy matrix 
discussed earlier. Longitudinal tensile tests were conducted until failure on composite laminates 
with various fiber orientations. Tests were conducted at strain rates of about 5xl0' 5 /sec, about 
1.0 /sec and about 400-600 /sec. Dog-bone shaped specimens were used with a gage length of 
approximately 0.9525 cm. The low strain rate testing was conducted using an Instron hydraulic 
testing machine. The high strain rate tests were conducted using a tensile split Hopkinson bar 
apparatus. More details on the experimental procedures can be found in Gilat, et al (ref. 19). 

The IM7/977-2 composite has a fiber volume ratio of 0.60. The material properties used in 
this study for the IM7 fiber include a longitudinal modulus of 276 GPa, a transverse modulus of 
13.8 GPa, a longitudinal Poisson’s ratio of 0.25, a transverse Poisson’s ratio of 0.25 and an in- 
plane shear modulus of 20.0 GPa. These properties are as given in Gates, et al (ref. 20), with the 
exception of the value for the transverse Poisson’s ratio which was taken from Murthy, et al (ref. 
21) based on representative carbon fiber data. The material properties of the 977-2 resin are 
listed in Table 1. In the computations, five fiber slices were used in the one-quarter of the unit 
cell which was analyzed. This value was found to yield sufficiently converged answers. 

Experimental and computed longitudinal tensile stress-strain curves for three laminate 
configurations ([10°], [45°] and [90°]) of the IM7/977-2 material are shown in Figure 7, Figure 8 
and Figure 9. In Figure 7, results for the [10°] laminate at strain rates of 1 .8xl0" 5 /sec, 0.56 /sec 
and 320 /sec are shown. In Figure 8, results for the [45°] laminate at strain rates of 4.75xl0' 5 
/sec, 1.2 /sec and 405 /sec are shown. In Figure 9, results for the [90°] laminate at strain rates of 
4.44x1 0' 5 /sec, 1.09 /sec and 405 /sec are shown. Note that the [10°] and [45 °] laminates display 
a strong rate dependency and nonlinearity in the deformation response across the entire range of 
strain rates. The experimental results for the [90°] laminates at low and medium strain rates 
show minimal rate dependence and nonlinearity. This lack of rate dependence in the [90°] 
laminates at lower rates is most likely due to the specimens failing abruptly due to phenomena 
such as fiber/matrix debonding while the effective stresses in the matrix were still primarily in 
the elastic range, and the matrix elastic modulus is the same at the low and medium strain rates. 
The fact that the computed results also predict a mostly linear defonnation response seems to 
support this hypothesis. At high strain rates, the stiffness and maximum stress of the [90°] 
laminate increase over the values observed at low and medium strain rates. This increase is most 
likely due to the increase in the matrix modulus observed at high strain rates. The apparent 
nonlinearity in the experimental stress-strain curve of the [90°] at high strain rates is most likely 
due to the lack of stress equilibrium at low strains and the specimen geometry. 

As can be seen in the figures, the analytical model captures the strain rate dependence (or 
lack thereof for the [90°] laminates at low strain rates) and nonlinearity observed in the 
experimental stress-strain curves. Quantitatively, overall the comparison between the 
experimental and analytical results is quite good. The stresses for the [45°] and [90°] laminates 
at high strain rates are somewhat under predicted at low strains. There are several possible 
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causes for this discrepancy. First, particularly for the [45°] laminate, recall that due to the nature 
of the tests on the bulk matrix the initial matrix elastic modulus at high strain rates was 
approximated. The modulus which was determined for high strain rate conditions might be too 
low. There also could be a lack of stress equilibrium in the experiments at low strains, which 
could lead to the experimental stresses being too high. The lack of stress equilibrium at low 
strains in the experimental tests is probably more significant for the [90°] laminate. This 
particular laminate orientation is dominated by the matrix. Recall that in the high strain rate 
tensile tests for this matrix material, in which a similar specimen geometry was used, the tensile 
stress-strain curves displayed a significant oscillation. If an improved specimen geometry was 
employed for the high strain rate tensile tests of the [90°] laminate, similar to the specimen 
geometry used in the high strain rate tensile tests for the PR520 resin, the oscillations in the 
experimental stress-strain curves would most likely be significantly reduced, and the match 
between the experimental and computed results would most likely be much improved. However, 
at higher strains for both the [45°] and [90°] laminates the experimental and computed results 
match reasonably well. In the prediction of composite damage and failure, which will be a future 
goal of this research, correctly predicting the composite stresses at higher strain levels will be 
important. 

To further explore the significance of properly accounting for the hydrostatic stress effects in 
the composite analysis, the tensile stress-strain curve for the [45°] laminate at the low strain rate 
of 4.75x1 O' 5 /sec is once again considered. In Figure 10, the experimental stress-strain curve, 
along with the original computed curve (“Alpha Included” in the figure), is presented. A tensile 
curve computed in which the hydrostatic stress effects are neglected (a o =ai=0, “Alpha=0” in the 
figure) is also shown. As can be seen in the figure, when the hydrostatic stress effects are 
neglected in computing the tensile response, the stresses are significantly over predicted as 
compared to the experimental results, similar to what was observed in the bulk polymer. These 
results indicate that the hydrostatic stress effects are significant for the analysis of polymer 
matrix composites as well as bulk polymers and must be accounted for within an analysis. 

Conclusions 

The Bodner-Partom viscoplastic state variable model (ref. 12) has been modified in order to 
analyze the strain rate dependent, nonlinear deformation of polymeric matrix materials. The 
effects of hydrostatic stresses on the inelastic defonnation have been accounted for by modifying 
the effective stress and effective inelastic strain definitions given in the original equations. The 
material constants are determined solely through the use of simple tension and shear tests. 
Furthermore, the effects of the hydrostatic stresses are assumed to vary over the course of the 
loading process, which leads to an improved correlation with experimental results. The tensile 
and shear defonnation response of two representative bulk polymers over a range of strain rates 
have been successfully predicted. The constitutive equations have been implemented within a 
strength of materials based micromechanics approach in which the unit cell is subdivided into a 
set of independently analyzed slices. The micromechanics technique was used to predict the 
strain rate dependent deformation of a representative polymer matrix composite for three fiber 
orientations. Qualitatively, the strain rate dependence and nonlinearity observed in experimental 
stress-strain curves was captured. Quantitatively, at low and medium strain rates the 
experimental stress-strain curves were accurately predicted by the analytical model. At high 
strain rates, there were some discrepancies at lower strain levels between the experimental and 
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computed results, most likely primarily due to issues in the experimental tests. Overall, the 
analysis methods presented here appear to correctly capture the important features of the 
nonlinear, strain rate dependent deformation response of polymers and polymer matrix 
composites. 

The analytical methods described in this study can be used in explicit finite element codes to 
provide a more realistic analysis of deformation during impact events for composite materials. 
The analytical model is computationally efficient so that it is feasible to use in the explicit finite 
element analysis of large composite structures. 
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TABLE 1.— MATERIAL PROPERTIES FOR POLYMER MATRIX MATERIALS 



Strain 

Rate/sec 

Modulus 

GPa 

Poisson’s 

Ratio 

D„ 

1/sec 

n 

Zo 

MPa 

Zi 

MPa 

q 

a G 

Oti 

PR520 

7x10° 

3.54 

0.38 

1x10 s 

0.93 

396.09 

753.82 

279.26 

0.568 

0.126 

1.76 

3.54 

420 

7.18 

977-2 

9x1 O' 5 

3.52 

0.40 

1x10 s 

0.85 

259.50 

1131.4 

150.50 

0.129 

0.152 

1.9 

3.52 

500 

6.33 



♦ Experiment-Low Rate 
“ Computed-Low Rate 
■ Experiment-Medium Rate 

— — Computed-Medium Rate 
▲ Experiment-High Rate 

- - - Computed-High Rate 


Figure 1. — Experimental and computed shear stress-shear strain curves for PR520 resin at strain 
rates of 7x1 O' 5 /sec (Low Rate), 1.76 /sec (Medium Rate), and 420 /sec (High Rate). 
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Strain 

Figure 2. — Experimental and computed tensile stress-strain curves for PR520 resin at strain rates 
of 5x1 O' 5 /sec (Low Rate), 1.4 /sec (Medium Rate) and 510 /sec (High Rate). 



Strain 

Figure 3. — Experimental and computed shear stress-shear strain curves for 977-2 resin at strain 
rates of 9x1 O' 5 /sec (Low Rate), 1.91 /sec (Medium Rate) and 518 /sec (High Rate). 
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♦ Experiment-Low Rate 
" Computed-Low Rate 
■ Experiment-Medium Rate 

— — Computed-Medium Rate 
▲ Experiment-High Rate 

- - - Computed-High Rate 


Figure 4. — Experimental and computed tensile stress-strain curves for 977-2 resin at strain rates 
of 5. 7x1 O' 5 /sec (Low Rate), 1.31 /sec (Medium Rate) and 365 /sec (High Rate). 



Strain 


Figure 5. — Effect of hydrostatic stress effect state variable a on tensile stress-strain curve for 

PR520 resin at strain rate of 5x1 O' 5 /sec. 
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Sliced Unit Cell 
Circular Fiber 



Sliced Unit Cell 
Rectangular Fiber Slices 



Slice of Unit Cell 


Portion of Sliced Unit Cell Analyzed 


Figure 6. — Schematic showing relationship between unit cell and slices for micromechanics. 



♦ Experiment-Low Rate 
" Predicted-Low Rate 
■ Experiment-Medium Rate 

— — Predicted-Medium Rate 
▲ Experiment-High Rate 

- - - Predicted-High Rate 


Figure 7. — Experimental and predicted stress-strain curves for IM7/977-2 [10°] composite at 
strain rates of 1.8xl0' 5 /sec (Low Rate), 0.56 /sec (Medium Rate) and 320 /sec (High Rate). 
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Strain 


Figure 8. — Experimental and computed stress-strain curves for IM7/977-2 [45°] composite at 
strain rates of 4.75xl0' 5 /sec (Low Rate), 1.2 /sec (Medium Rate) and 405 /sec (High Rate). 



♦ Experiment-Low Rate 
- Predicted-Low Rate 
■ Experiment-Medium Rate 

— — Predicted-Medium Rate 
a Experiment-High Rate 

- - - Predicted-High Rate 


Figure 9. — Experimental and computed stress-strain curves for IM7/977-2 [90°] composite at 
strain rates of 4.44xl0' 5 /sec (Low Rate), 1.09 /sec (Medium Rate) and 405 /sec (High Rate). 
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Strain 


Figure 10. — Effect of hydrostatic stress effect state variable a on tensile stress-strain 
curve for IM7/977-2 [45°] laminate at strain rate of 4.75xl0' 5 /sec. 
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